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O ■ Abstract 

, This paper investigates the ability of the stochastic subspace identification technique to return a 

CO ■ 

■ valid model from finite measurement data, its asymptotic properties as the data set becomes large, and 

asymptotic error bounds of the identified model (in terms of H2 and Hoo norms). First, a new and 
straightforward LMI-based approach is proposed, which returns a valid identified model even in cases 
where the system poles are very close to unit circle and there is insufficient data to accurately estimate 
the covariance matrices. The approach, which is demonstrated by numerical examples, provides an 
altenative to other techniques which often fail under these circumstances. Then, an explicit expression 
for the variance of the asymptotically normally distributed sample output covariance matrices and block- 
Hankel matrix are derived. From this result, together with perturbation techniques, error bounds for the 
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■ state-space matrices in the innovations model are derived, for a given confidence level. This result is 

in turn used to derive several error bounds for the identified transfer functions, for a given confidence 
level. One is an explicit bound. Additionally, two error bounds are derived; one via perturbation 
analysis, and the other via an LMI-based technique. 

«- , 

X 

I— 1 . Index Terms 

a 

Asymptotic variance, stochastic subspace identification, positive realness, H2 norm error bound, 
Hoo norm error bound, linear matrix inequalities. 

Q. Li is currently with the Department of Electrical Engineering and Computer Science, Massachusetts Institute of Technology, 
Cambridge, MA 02139 USA (e-mail: quanli@mit.edu). 

J.T. Scruggs is with the Department of Civil & Environmental Engineering, University of Michigan, Ann Arbor, MI 48109, 
USA (e-mail: jscruggs@umich.edu). 



December 27. 2012 



DRAFT 



2 

I. INTRODUCTION 

In system estimation and control design, an undoubtedly crucial issue is the identification 
of the stochastic part of the plant JT], [0, [|3). To obtain a nominal model for the stochastic 
system, stochastic subspace identification methods are often used flU, 0, J6l. However, due to 
model error or data insufficiency, these methods may encounter a failure mode without any valid 
model returned. This is especially common when the system poles are very close to the unit 
circle, as demonstrated experimentally in [7j. To overcome this difficulty, several approximation 
approaches for guaranteeing a valid model have been proposed in recent years flU, [0, [fTOl . 
ifTTfl. In jH, Mari et al. proposed an improved version of stochastic subspace identification 
algorithm in which linear matrix inequality (LMI) based techniques are used to constrain the 
identified system poles inside the unit circle as well as for multivariate covariance fitting, thus 
guaranteeing the solvability of the associated discrete algebraic Ricatti equation (DARE) and 
a valid model returned. Although promising, this method relies on coprime factorizaiton for 
covariance estimation. For large-dimensional multiple input multiple output (MIMO) systems, 
the numerical robustness of coprime factorizations may be problematic, and in any case is not 
well understood. Goethals et al. in [|9) used regularization techniques to impose positive realness 
on the associated covariance model, which we will define in Section 11(A). The solvability of 
the associated DARE and thus the feasibility of a valid model are then satisfied, although at 
the cost of introducing a small bias on the identified model. Different from these approaches, 
the present paper first exploits an equivalence between the solvability of the DARE and the 
nonemptyness of a convex set, and then by the Positive Real Lemma, establishes in section 
11(B) a more straightforward approximation approach based on LMI techniques. 

To estimate how much data is required to identify a model given a model error bound 
(such as "Hoc, norm bound) with a chosen confidence level as a starting point for further robust 
controller design, the asymptotic analysis and perturbation methods provide fundamental tools. 
For identification of the deterministic part of the plant, asymptotic statistical properties for 
prediction error (PE) methods |{I2]|. instrumental variable (IV) methods lfT3l . and subspace 
methods lfl4l . are all well-established. However, the corresponding asymptotitc properties are 
not thoroughly studied for the identification of the stochastic part of the plant. 

Due to a finite sample size and influence of system and measurement noise, a model error 
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always exists between the identified and true system. To assess the quality of the identified 
model, a model error bound in terms of H 2 or Hoo norms is often useful. Moreover, robust 
control theory is often predicated on knowledge of an upper bound for the norm of the 
model error. To derive an norm model error bound for the deterministic part of a system, 
various identification methods are proposed by [[151 . [fl6l . [fTTl . [fT8l . |fl9l . In |fl9ll , an upper 
bound on the "Hoc norm of the model error is derived via a frequency response curve fitting 
procedure which minimizes a maximum amplitude criterion and guarantees the stability of the 
identified model. In this procedure, linear as well as nonlinear programming techniques are 
used. Another approach to quantify a model error, proposed in ll20ll . establishes a framework 
connecting PE methods with robust control theory. In this framework, the tools of PE methods 
are used to quantify an uncertainty region to which robustness tools are conveniently adapted 
such that robustness analysis of a controller and the quality assessment of the uncertainty region 
are easily carried out. For the stochastic part of a system, the quantification of model error is 
rarely reported. One of the contributions of this paper is to derive a model error bounds in terms 
of both H 2 and Hoo norms, with a confidence level given by asymptotic analysis of stochastic 
subspace identification. 

In this paper, we first propose a new and straightforward LMI-based optimization approach 
to identifying the stochastic system. Then for the identified stochastic system, we investigate its 
asymptotic behavior using asymptotic analysis and perturbation methods. Particularly, we derive 
the asymptotic Frobenius norm (F-norm) error bounds of the state-space matrices in innovations 
model. With these asymptotic properties, we derive the explicit expressions of the H 2 an d %oo 
norm model error bounds for the identified system, for a given confidence level. We also propose 
an LMI-based approach to computing the norm model error bound. In order to verify our 
analytical results, we first apply them to identify an innovations model for a plant with poles 
are very close to unit circle, and for which the amount of measurement data is insufficient to 
accurately evaluate output covariances. Then, we propagate these derived bounds, to estimate 
H 2 and "Hqo model error bounds for the identified system. Based on simulation results, it is 
shown that the improved stochastic subspace identification procedure guarantees a valid model 
returned, together with 1-L 2 and norm model error bounds with a chosen high confidence 
level. 

The outline of this paper is as follows. In Section II, the stochastic subspace identification 
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procedure will be outlined, and the equivalence of its failure mode with the emptyness of a convex 
set will be illustrated, which motivates an LMI-based optimization approach for guaranteeing a 
valid model returned. The asymptotic analysis and perturbation method will be given in Section 
III, deriving the asymptotic distributions of state-space matrices in the associated covariance 
model, and the asymptotic F-norm error bounds of the state- space matrices in innnovations 
model with a chosen confidence level. In Section IV, these asymptotic results are combined to 
derive the I-L2 and Hoo norm model error bounds of the identified system with a chosen confidence 
level. Two typical numerical examples will be presented in Section V, and the conclusions will 
be drawn in Section VI. For brevity, hereafter we use / or I n E TZ nxn to denote the identity 
matrix with compatible dimension, (■) to denote the estimated value of (•) from the measurement 
data, 6(-) = (■) — (•) to denote the perturbation of (■) due to finite data samples, o(-) to denote 
the higher order perturbation satisfying o(-)/(-) — > as (•) — > 0, and = to denote the first-order 
approximation due to the perturbations from the finite sample, vec(-), <g>, and tr(-) represent the 
columnwise vectorization of a matrix, Kronecker product and the trace of a matrix, respectively. 
(•)*, (-) T and (-) H represent conjugate, transpose and conjugate transpose, respectively. || • \\p 
and || ■ || 2 represent F-norm and 2-norm, respectively. 

II. AN IMPROVED STOCHASTIC SUBSPACE IDENTIFICATION ALGORITHM 

In this section, we will first recall the standard stochastic identification flU, [[HI in which 
the associated DARE may be unsolvable especially when measurement data is insufficient to 
estimate stationary output covariances, and the system poles are very closed to the unit circle. 
To overcome this difficulty, we establish an equivalence between the solvability of the DARE 
and the positive realness of an associated covariance model, and then propose an LMI-based 
approach to impose positive realness on the covariance model, thus guaranteeing a valid model 
is returned. 

A. Standard subspace identification procedure 

Assume the following state- space innovations model is a minimal realization of the vector 
stochastic process {y k } 



Xk+l — 



Ax k + KQ 2 e k 



y k = Cx k + Q 2 e k 



(1) 



December 27, 2012 



DRAFT 



5 



where %k G lZ nxXl ; yk G 1Z nyXl ; K is the Kalman gain; Q G TZ nyXny is the innovations covariance 
matrix; G 7£. n!/><1 are the nomalized white innovations with covariance matrix E[eke1] equal 
to I n . We assume zero-mean processes throughout. For k G Z> , denote the process covariances 
as 



R k = E\y(t)y 1 (t - k)} 



(2) 



According to the minimal realization assumption, the following observability matrix and 
controllability matrix T are both in full rank. 

C 



n 



CA 



CA 



m—l 



T=[D AD ■■■ A m - l D\ 



(3) 



where D = E[xk+\y1] an ^ m > n x- Noting that the sequence {R{\ are the Markov parameters 
of the system {A, D, C, Ro}, the following factorization of the block-Hankel matrix of the output 
covariance holds: 

Ri i?2 • ■ ■ R m 

R2 R3 ' ' ' Rm+l 



H 



nr 



(4) 



Rm Rm+l R2m-l 

From © and ©, the rank properties of and T imply that 

rank(if) = n x , for m > n x (5) 

We refer to (A, D, C, Rq/2) as the covariance model associated with the innovations model ©. 
From singular value decomposition (SVD) of H, we have 

A s 




nr = [u s u n ] 



[v 3 v n ] 



u s k s v s 



(6) 



where the diagonal matrix A s G R nxXnx is nonsingular. © indicates the realizations of f2 and 
r as 

n = u s aIt, r = t^aIv* (7) 
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where T G TZ nxXrix is any (nonsingular) similarity transformation matrix. Hereafter we set T = I. 
In view of ©, it is indicated that C is the first n y rows of f2, and D is the first n y columns of 
T. A can then be solved by the following overdetermined linear equations: 

VLA = (8) 

where f2 and O are respectively the first m — 1 matrix blocks and the last m — 1 matrix blocks 
of fi. We then solve the discrete algebraic Ricatti equation (DARE) to obtain a positive definite 
covariance matrix P = E[xkxj?\ 

P = APA T + (D - APC T )(R - CPC T Y\D - APC T ) T (9) 

and K and Q are then given by 

Q = R - CPC T 

K = (D — APCf r )Q~ 1 (10) 

B. Correcting for errors due to finite measurement data 

In practice, R k is estimated by the empirical (i.e., statistical) sample output covariance 

1 N 

t=k+l 

The above algorithm is then implemented assuming R k The resultant estimations for 

{^4,(7,(5,^} are denoted {A, C, Q, K}. As N — > oo, the estimations converge to the true 
values, but for finite, the inevitable existence of estimation error in R k leads to two problems. 

The first potential problem is that A, as determined as the least-squares solution to d8}, may 
be unstable. If this is the case, a simple LMI-based correction proposed in flU can be used to 
recover stability: 

min Aw \\(A-A)W\\ F 

s.t. W - AWA T > 

W > (12) 

where A is the adjusted (and asymptotically stable) approximation of A. 

The second potential problem is that the DARE © may fail to yield a positive-definite, real 
solution P, thus prohibiting the derivation of approximate innovations model parameters {Q, K} 
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as in (fTQb ■ Failures are especially likely for large dimensional systems with poles close to the 
unit circle. The following result reveals the equivalence between the solvability of DARE and 
the positive realness of the covariance model. For the theorem, we assume a covariance model 
{A, D, C, -Ro/2} has been derived via the procedure above. 

Proposition 1: Given a symmetric, positive definite matrix R and a controllable and ob- 
servable covariance model {A, D, C, R /2}, the associated Riccati equation © has a positive 
definite solution P with innovations covariance Q = Ro — CPC T > if only if the set V, 
defined as 

P-APA T D-APC T 
D T -CPA T R n -CPC T 



V 



P\ 



> 0,P > 



(13) 



is nonempty. 

Proof: flU and llTOl have shown that given a controllable and observable model, © has a 
positive definite solution P with R — CPC T > if and only if C(zl — A)~ 1 D + R /2 is positive 
real. Thus, we have D T (zI — A T ) _1 C* T + Ro/2 is also positive real. Applying the Positive Real 
Lemma obtains the conclusion. ■ 
Thus, the failure of DARE indicates a null set V resulting from the estimate (A,D,C,Rq). 
This suggests that a feasible member in V be approximated by solving the following optimization 
problem 

P-APA T I) M'C <I>,, <I> 

D T -CPA T R n -CPC T 



$11 $12 
$22 



^2 



arg mm 

P,«I>ll,<I>12, < I > 22 



12 



$1 2 $22 



S.t 



$H $i 2 
$f 2 $22 



>0 



P >0 



(14) 



in which we can use 2-norm or F-norm to express it as a semidefinite program. For simplicity, 
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choosing 2-norm in ([141) gives the following LMI problem 



min A 

AI P-APA T -$ U D-APC T -$ 12 

AI D T -CPA T -^ 2 R -CPC t -^ 22 

I 

(sym) I 



> 



(15) 



P > 



$12 
$22 



> 



(16) 
(17) 



where the variables are the positive definite matrix P, the nonnegative definite matrices $n and 
$22, the non-symmetric matrix $i 2 , and the scalar A. 

For convenience, we keep the estimated A and C unchanged, and adjust D and R for 
guaranteeing that a positive definite solution to © exists. After solving (fl4l . we obtain P 
as the solution to the following discrete Lyapnov equation: 



P = APA T + $n 



(18) 



and then the new D and R are adjusted from D and R by 

D = APC T + $i 2 (19) 

R = CPC T + $ 22 (20) 

Substituting the adjusted D and R to DARE ©, it is guaranteed that DARE has a positive 
definite solution P, and finally, the estimated K and Q are updated by 

Q = R - CPC T (21) 

K = (D — APC^Q- 1 (22) 

We thus arrive at adjusted parameters {A, K, C, Q} as the approximation to the innovations 
model parameters from (Q3. 
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III. ASYMPTOTIC AND PERTURBATION ANALYSIS OF STOCHASTIC SUBSPACE 

IDENTIFICATION 

In this section, we will derive the asymptotic distributions of the empirical sample output 
covariance matrices, with its variance expressed in terms of the power spectral density (PSD) 
of the vector stochastic process y^. Later we will use matrix perturbation analysis to derive 
the asymptotic variances of the asymptotically normally distributed state- space matrices in the 
covariance model, and the perturbed state-space matrices in the innovations model in terms of the 
perturbations 8Rq and 5H . Using Chi-squre cumulative distribution function, the F-norm error 
bounds of state-space matrices in the innovations model are given with a chosen confidence 
level. 

A. Asymptotic Distributions of the Sample Output Covariance Matrices 

To facilitate the deduction of asymptotic distributions of the empirical sample output covari- 
ance estimates of the vector stochastic process {yk}, we introduce a Lemma which illustrates 
that the expectation of the mixed product of scalar and vector random variables can be expressed 
in terms of first- and second-order moments. This Lemma, as a special case of theorem 1 in 
ETi . gives a useful tool to calculate the fourth order moment. 

Lemma 1: If Xi,x 2 G 1Z and X 3 ,X 4 G 1Z nxl have jointly gaussian distributions, then 



We can now propose a theorem which reveals the asymptotic distributions of the empirical 
sample output covariance estimates of the vector stochastic process {yk}- 

Theorem 1: Consider the empirical sample covariance .R for the vector stochastic process 
{Vk} given by (OQ) and the empirical Hankel matrix H given by substituting each R k in © by 
R k in (fTT|) . Then R and H are both asymptotically normally distributed 



E [ Xl x 2 X 3 Xj] =E [x x x 2 ] E [X 3 X 4 T ] + E [ Xl X 3 ] E [x 2 Xj] + E [x 2 X 3 ] E [ Xl Xj] 
- 2E[x 1 }E[x 2 ]E[X 3 }E[Xj] 



(23) 




(24) 
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where the variance matrix Vrq,h 
P H given by 

1 



Pr Pr h 
p RqH p H 



with the block matrices Pr , PrqH and 



Pro ={Inl + K ny )— / S* y (u) ® S , ,(o;)da; 



2tt 
1 

2tt" 



<8> ^(w))* g> (Ef g) 5 y (w)) do; 



(25) 
(26) 



1 



=(4> 2 + ^n,™)— / ((E 2 Ef ) ® ® ((ExEf) ® do; (27) 

where S^w) is the PSD of the vector stochastic process {y k }; K n is a permutation matrix 
satisfying 



K n (e { <g> ej) 






(28) 



; Ex = [1 e iw ... e ^ m -V] T ; E 2 = [e ~ < — ... < 



, is a n-dimensional unit vector with 



the zth element be 1 and others 0. 

Proof: We first proof (|23T ). and then extend the proof to obtain (12~6l and (1271) . We have the 
covariance of vec(_R — Ro) 

Cov(vec(i? — Ro), vec(i? — .Ro)) 

= ^ J] ^ [y k ® y k - vec{R )) (y h ® y h - vec(R )) T (29) 

l<fc,h<iV 

where we derive the second line using the fact vec(y k yl) = y k ®y k . Each entry in the sum in 
d29b is 



E 



(yk <&yk~ vec( J R )) (Vh ®Vh~ vec(i? )) J 



-E [(y k yl) ® [Vkvl)] - vec(i? )vec( J R ) T 

= ^ej ® £ (yfy^yuyl) - vec(R )vec(R 



+ X> e i ® E(y%\ k )E(yfy£) - vec(i? )vec(i? ) T 



(30) 
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where we use Lemma \T\ to derive the fourth line from the third line. For the first term in (|30l ), 
we have 



hi 

= (e ? ej) ® (e m e T n ) E^y^E^y^) 

i,j,m,n 

= E e ™) ® e «) EtffyPmyMyM) (31) 



i,j,m,n 

For the second term in (|30T ). we have 



Y^®^(yfy k )E(y^yl)) 

y 

= E( e '® ( e i ® EivPvh))] 

= Y-<E{yfy k )e^ C \E{y^y h ))e T J ) 

= vec(E(y k yl))vec T (E(y h yj;)) 

= vec(R )vec(R ) T (32) 
For the third term in (|30l ). we have 

i,j,m,n 

E ((e m ej) ® (*«£)) E{$ y ®)E( y ^) 

i,j,m,n 

E (e™ ® e„) (ej ® e T n ) E{yfy^)E{y^) 

i,j,m,n 

K ny E (ei ® e m ) (ej ® e£) E{y^yf)E{y^) (33) 



i,j,m,n 
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where we switch the indices i and m to derive the third line from the second line. Substituting 
(HE ~ 03b to we have 



E 



(yk ®Vk- vec(i? )) {Vh ®Vh- vec(Ro)) 



= (7„2 + if n j£ ® £ (y k yl) 

= (J n 2 + K ny )R k „ h <S> Rk-h 
Substituting (O to (l29l> . we have 

Cov(vec(i? — -Ro)? vec(i? — i? )) 

l<k,h<N 

= l(I n2x+ K ny )J2R(r)®R(r) 

T 

= hl n% + K ny )±- r J2 ® ^(r 2 )e--(— du; 



(34) 



do; 



= — (l n 2 + K n ) — 



where S y (u;) denotes the PSD of the vector stochastic process {yk}- Thus, (|25l) is concluded. 
For P ff , following the same procedure we have 

Cov (vec(H - H),vec(H - H) 



( 1 N 1 N \ 

Cov vec(- Vi,kV 2 T k - H),vec(— £ V 1>k V* k - H) 

\ k=l k=l J 



Sy (u) g) S Vl (u)duj 



(35) 



where V x , k = [y k y k+1 - vl+m-if and V 2,k = [vLi vl-2 - vl-J T \ 5 VxH and S V2 (u) denotes 
the spectral densities of V\ k and V 2k respectively, and can be obtained by 



S Vl (u) = (E 1 ® In.) S y (u) (Ef ® I nx ) 
= (E X E*) ® S v {u) 



(36) 
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Likewise, we have 

S V2 (uj) = (E 2 <g> J n J S y (u>) (E* ® J n J 

= {E2E2) ® S y (u) (37) 

Substituting (1361) and (1371) to (1351) gives (1271) . For PrqH, following the same procedure we have 
that 



Pr h =Cov (vec(i? - Ro), vec(# - #)' 

/ N \ 

=Cov vec(- £ y fcJ £ - i2o), vec(- £ V l>k V? k - if) 

V fc=l k=l J 

=^(I nl + K ny )±- j S* yV2 (u) ® S yVl (u)du (38) 

where S y y 2 = E2®S y (uj) is the cross spectral density of {y k } and {V^fc}, and S y v 1 = E^^S y (u) 
is the cross spectral density of {y k } and {Vi,*}. Likewise, we have ((2~6"1) . Proposition 7.3.2 ~ 
7.3.4 in |[22~1|. as well as theorem 1 in [23], claim that the sample covariances R k in (fTT|) are 
normally distributed as the measurement data size N goes to infinity. Thus, (124)) is conluded. ■ 
Remark 1: the PSD S y (u) of vector stochastic process {y k } can be computed by the transfer 
function G e (z) from the normalized innovations to output, i.e. 

S v (u) = G e (e iuj )G T e (e' iuJ ) (39) 

where 

G e (z) = C(zl - Ay 1 KQ^ + (40) 



B. Asymptotic Distributions of the State-space Matrices in the Covariance Model 

Based on the asymptotic distributions of the empirical sample covariances of the vector 
stochastic process {y k }, perturbation analysis of SVD is firstly applied to derive the asymptotic 
distributions of the controllability matrix T and observability matrix f2, from which we then 
derive the asymptotic distributions of the state space matrices in the covariance model. 

Assume the true covariance matrix H has SVD 

H = U S A S V S T + U n A n V^ (41) 
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where the diagonal matrix A s G TZ nxXnx is defined in © and A„ = 0. Due to finite data samples, 
there exists a perturbation 5H in H 

H = H + 5H (42) 

such that there exist the corresponding perturbations in the subspaces and singular values. The 
SVD on H gives 



H = U S A S V S T + UnAnV? 



(43) 



Due to 5H, all the terms on the right hand side of (|43l may differ from those on the right hand 
side of (|41"I) . [24J developed the perturbation analysis of SVD to the second order. We will apply 
its main theorem to analyze how 5H influences the SVD of H. Accordingly, we have that 

A s = A s + Uj5HV s 

U S = U S + U n Ul5HV s K~ s l 

V S = V S + V n V^5H T U s A; 1 (44) 
Also, matrix perturbation analysis obtains the first order approximation of AJ 

A! = AJ + 5A sq (45) 

where 



vec(<5A 



A! +AJ 



(V S T ® Uj)vec(8H) 



(46) 



Algebraic manipulation yields the perturbations of the controllability and observability matrices 
in the covariance model in terms of 5H 



vec ( UAl - UM 



s 1 *-s 

_ 1 _ 



vec [AIV S J -AIV S ' 



nivec(5i?) 
= U 2 vec(5H) 



(47) 
(48) 



where 

n 9 



{In y ® U s 
(K ® In y . 



A| +AJ 



/_ ® A! + A! 



5 (£4^ T ) 

a:% t 



(49) 



From theorem [Hand (l47l ~ (|48|) . we have the following proposition which reveals the asymptotic 
distributions of the observability and controllability matrices in the covariance model. 
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Proposition 2: Consider the estimated observaility matrix Vt and controllability matrix f given 
by © for T — I nx . They are asymptotically normally distributed 



'Nvec{Q -ft) — > Af(0, Pq) 
v / iVvec(f - T) — > AT(0, P r ) (50) 
The asymptotic variance matrices are 

P n = U.PnUl (51) 
P r = n 2 P H II^ (52) 

where Hi and H 2 are given by (|49i 

By this proposition, we have the following theorem for the asymptotical distributions of the 
state-space matrices in the covariance model. 

Theorem 2: Assume m > n. The state-space matrices A, C, and D estimated from Section 
11(A) are asymptotically normally distributed 



where 



'Nvec{A -A) — y Af{0, P A ) 
\//Vvec(<5 - C) — )• jV(0, P c ) 

VNvec(D -D) — > Af(0, P D ) (53) 
P A = EP n E T , P c = (I nx <g) $ 3 )Pn(/^ ® $ 3 ) T 

P D = ($J <g> /njP r ($I ® /nJ T (54) 



and 



$ x =[/ ; 0] G Jl( m - 1 ) n y xmn y ( $ 2 — [0, J] G Jl( m - l ) n v xmn y 
$ 3 =[J > 0] G K nyXmrH >, $ 4 = [J, 0] T G ^ mn ^ xn ^ 

s=/ n;c ® ^(n T $f$ 1 n) _1 n r $f$ 2 ) -a t ® ((n T $f$ 1 n)" 1 j2 T $f$ 1 ) (55) 

Proof: From the structure of ®, we have 

vec(<5 - C) = vec($ 3 5fi) = (7 n , ® $ 3 )vec(<m) 

vec(P - £>) = vec(5r$ 4 ) = ® I n Jvec(6r) (56) 
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From ([8]), the least-square solution of the overdetermined linear equation for A is 



-i 



a= ($iny r $ 2 ft (57) 

First-order approximation of the right-hand side in (1571 and then vectorizing A — A yields 

vec(i - A) = EvecSQ (58) 
From (|56]) and (l58i applying © yileds (l53l . ■ 



C. Perturbation Analysis of the State-space Matrices in the Innovations Model 

To facilitate the perturbation analysis, the general DARE © is transformed to its equivalent 
form. 

P - A T S P(I - SP)- l A s - M = (59) 

where M = DR 1 D T , A s = A T — C T R l D T and S = C T R 1 C. In terms of 5A S , 5M and 5S, 
G51 gives the first-order perturbation SP 

vec(5P) = J^vec(5M) + Jf 1 J 2 vec(5A s ) + J^J 3 vec(5S) (60) 

where Ji = I nl - A T C ® A T C , J 2 = (A T C P <g> I nx )K nx + I nx ® A T C P, J 3 = A T C P g> 4fP and 
A c = (I — SP) -1 ^. Also from the expressions of A a , M and S, we have that 

vec(5M) =J 4 vec(5P ) + </ 5 vec(5P) 

vec(54,) =lf na; vec(<L4) + J 6 vec(5C) + J 7 vec(5R ) + J 8 vec(<5L>) 
vec(55) =J 9 vec(<5P ) + J 10 vec(5C) (61) 

where 

J 4 = -{DRo 1 <8> DRo 1 ), J 5 = {I nx ® DR Q v )K n ^ ny + P^ 1 <g> /„, 
J 6 = -(PPq- 1 <8» InjK^, J 7 = PP 1 ® C r Po \ = -(4, ® C T R^)K n ^ y 
J 9 = -(C^ 1 ® C^ 1 ), J 10 = (C T P - 1 ® I nx )K ny>nx + J n:c <g> C T P 1 (62) 
Substituting ([61) to d60]) obtains 

vec(5P) =J 1 " 1 (J 4 + J 2 J 7 + J 3 J 9 )vec(SR ) + J{ 1 J 2 K Ux vec(5A) + Jf 1 (J 2 J 6 + </ 3 </io)vec(<5C) 
+ Jr 1 (^5 + ^2 J 8 )vec(5P) (63) 
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where from (07), (g8), ([56) and (1581 . we have that 

vec(5A) =Hili[0 m 2 n 2 g J„ 



vec(5i?o) 
vec(5H) 



vec(SC) =(/ n;c ®<l>3)n 1 [0 



m 2 n 2 ,n 2 -^m 2 n 2 ] 



vec(5£>) =($J ® 4jn 2 [0 
Then, using (l64l) we express (|63l) as 

vec(5P) = 0! 

where 



m 2 ri 2 ,n 2 -^m 2 n 2 I 



vec((5i? ) 
vec(5H) 

vec(5Ro) 
vec(5H) 



(64) 



vec(5i?o) 
vec(5H) 



(65) 



-"m 2 ?! 2 ^ 2 -'m 2 n 2 J 



(66) 



£l — <A (<^4 + ^2^7 + <^3^9) [-^n 2 n 2 im 2„2] + JfVaATn.SnitO 

+ J\ {J2J6 + J?,Jl0)(In x <S> $ 3 )ril[0 m 2„2 in 2 J m 2, 
+ ^(^S + ^2^8)(^I ® 4jn 2 [0 m 2„2 in 2 I m 2„2] 

Denote B = KQ2 and F = Q^. Next, we will derive the explicit expressions of the perturbations 
5B and 5F. From (flOl ). applying perturbation analysis obtains 

vec(c5P) 

= (J®Q* +Q> ®I) _1 [vec(Mk) - (CP ® I + (I <Z> CP)K n ^ ny )vec(5C) - [C ® C)vec(SP)] 
vec((5P) 

= (0^5 ® J)vec(£D) - (Q'^CP ® J)vec(M) - (Q~*C ® A)vec(5P) 
- ® AP)2f n ., nw vec(*C) - (Q~* ® P)(/ ® + g, /)- 1 V ec(5Q) 
Also from Q = R — CPC T , applying perturbation analysis obtains 

vec(SQ) = vec(5R ) - (C ® C)vec(<5P) - [(J ® CP)K n ^ ny + {CP ® I)]vec(<5C) 

Substituting (|64h . (1651) and ([68]) to ([67]) yields 

vec(5i?o) 



(67) 



(68) 



vec(5F) =(£ 2 - <?3<?i ' 
vec(£i?) =(£ 4 + 



vec(<Jif) 

vec((5i? ) 
vec(5H) 



(69) 
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where 

g 2 = (i® g 1 / 2 + q 1 ' 2 ® j)" 1 {[j ng o ngxm2ni ] 

- (CP ® /+ (/ ® CP)^,^)^ ® $ 3 )n i [0 m 2„2 x „2 / m 2„2]} 

£3 = (/ ® Q 1/2 + Q 1/2 ® I) ~* (C (8) C) 

£4 = (Q~* ® *)($4 ® 4jn 2 [O m 2 n 2 x „2 / m 2 n 2] 

- (Q~lCP g) J)SIll[O m 2 n 2 xrl 2 J m 2„2] 

- g) AP)K nxiTly (I nx g> $ 3 )ni[0 m 2„2 x „2 J m 2„2] 

- (q-3 ig) ® g 1 / 2 + g 1 / 2 g, j)-!^ o n 2 Xm2ri 2] 

+ (g-i ® p)(j ® g 1 / 2 + g 1 / 2 ® i)- 1 ^ ® cp)K nx , ny 

+ (CP ® /)}(/„„ g) $3)ni[O m 2„2 x „2 J m 2„2] 

£5 HQ'** ®B)(I® Q 1 ' 2 + Q 1 ' 2 ® J)- 1 ^ <8) C) - (Q-5C <8) A) (70) 

According to theorem [T] as the data size N is sufficiently large, we can approximately quantify 
the F-norm error bound of the state-space matrix perturbations 5 A, 5B, 5C and 5F by Chi-square 
cumulative distribution function with a given confidence level. For SB, we have that 

\\5B\\ 2 F =vec(5B) T vec(5B) 

where g is a normally distributed vector with variance equal to I n 2 +m 2 n 2, i.e. g — > A/"(0, l n \+rrfln\)- 
Thus, we have 

g T g < Xl (72) 

with the confidence level given by Chi-square cumulative distribution function a(n 2 + m 2 n 2 , x 2 ) 
at each of the values in x 2 a using n 2 + m 2 n 2 degrees of freedom. Likewise, with this confidence 
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level, it is claimed that 



\\sb\\f 


V 2 

~ N 




) T (Qi + 






11*41! 


- N 


l^!ff(Sni[Om»n= 


xn 2 An 


2 n^]) T ( Sn i[O m 2„2 xn 2 




II^|If 


~ N 




III [O m 2 


n y xn y -^m 2 n 2 ]) i.{fn x @ 


S$ 3 )IIi[0 


11**11! 


~ N 


\v]l 2 H {g-2-g^i 


f(Q2 - 







1/2 



,1/2 



(73) 
(74) 
(75) 
(76) 



Remark: provided that © fails due to insufficient data, (fl4l) . as well as <fT9~l) and (|2"0l) . is required 
to adjust D and R for a valid model. In this case, we shall estimate the perturbations 5D and 
5R by the following equations 



5D = (D -£>) + (£>- D) 
5R = (R — Rq) + (Rq — Rq) 



(77) 



where D — D and Rq — Rq result from the LMI-based adjustments in (1191) ~ (1201) . and D — D 
and Rq — Rq are asymptotically normally distributed, as shown in (|53l and < f24b . 



IV. % 2 AND NORM MODEL ERROR BOUNDS 

In this section, we will derive the "H 2 an d "Hoo norm model error bound with a given confidence 
level. For brevity, we denote the true innovations model transfer function and the identified one 

as 



GJu) 



A 


B 


C 


F 



A 


B 


C 


F 



(78) 



A. H,2 Norm. Model Error Bound for the Identified System 

Theorem 3: Let (A,B,C,F) and {A,B,C ,F) be state-space representations of the original 
and identified systems, as shown in (1781) . such that 



(79) 



5A = o(A), 5B = o{B) 
5C = o{C), 5F = o{F) 



Assume, without loss of generality, that the state matrices A and A are Hurwitz, and that the 
data size N is sufficiently large such that (|24) approximately holds. Then, with the confidence 
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level a(riy + m 2 n 2 , x«), the V, 2 norm of the error system can be bounded by 

\\G e (u) - G e (u)\\l 2 <\\8F\\l + \\5B\\ F \\P\\ F + 2||B|| F ||<JB|| F ||<JP 1 || F + \\B\\ F \\5P 2 \\ F (80) 



where 5P\, 5P 2 , 5A, 5B, 5C and 5F are bounded by 



\\SP 2 \\ F <\\(A T ®A T - J 4n| )- 1 || 2 (2p|| F || ( 5P 1 || F ||M|| F + \\P\\ F \\5A\\l + \\5C\\ F ) 



||^||f ^-^ll^iio,if( Sn i[ m 2 n2xn2 On 2 ,]) (SIIi [0 m 2 n 2 xn 2 OngD^i^H II 2 



1/2 



iV 

II^IIf <§ I|p1h(& + G 5 Qif(G, + 060i)PjUl| 2 



(81) 
(82) 
(83) 
(84) 



■Xc 



\\5cf F <^\\v)H H {(i nx ® $ 3 )n 1 [o m2 „2 Xn 2 i m ^fX(i nx ® $ 3 )n 1 [o m2n 2 Xn? / manS ])^|| 2 

(85) 
(86) 

(87) 



1/2 



¥F\\l <^\\pLh(G 2 - Q,Qif{Q2 - GzGi)Pl, H h 

, P is the positive definite solution of the following Lyapunov equation 

A T PA -P + C T C = 
, and 





A 







B 


A = 






, B = 









A 




B 



, C=[C -C] 
M 1 = - [A T ® A T — P] 1 { [{A T P ® I)K 2nx + (J <g> A T P 
+ [(I ® C T ) + (C T ® I)K nyMx ] 

Proof: For brevity, we denote 



( 










) 






® 
























( 











v y y J 
















5 A = 



















55 



, <SC = [0 - <JC] 



(88) 



(89) 
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Consider the transfer function of the error system 





A 





B 


G e (u) - G e (u) = 





A 


B 




c 


-C 


F-F 



(90) 



The V.2 norm of G e (u) — G e (oj) is computed by an algebraic approach 



\GJu)-GJu) 



\%2 



tr 



(F-F) T (F-F)+ [ B J V J 



11^111 + tr 



- " - 

5 / B 



(91) 



(92) 



(93) 



where P is the solution to the following discrete Lyapunov equation 

(A + 5A) T V{A + 5 A) -V + (C + 5C) T {C + 5C) = 
Assume the asymptotic expansion of the solution V is 

V = P + 6P 1 + 5P 2 + o(6P 2 ) 

where 5 Pi = o(P) and 5P 2 = o(5Pi). Applying dominant balance to the perturbed equation 
(HU) yields 

A T PA -P + C t C = (94) 
A T 5PiA - 5 Pi + 5A T PA + A T P5A + C T 5C + 8C T C = (95) 
A T 5P 2 A - 5P 2 + 5A T 5PiA + 5A T P5A + A T 5P X 5A + 5C T 8C = 

From (|94| ), we have that 

r X -X 

-X X 

where X is the solution to the following discrete Lyapunov equation 

A T XA -X + C T C = (98) 



(96) 



(97) 
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Denote 8 Pi as 



8Pi 



8Xu 8X12 

^22 



8X12 8X-). 



Substituting d97J) and d99j to d95j) yields 

A T 5X n A - 5X n = 

A T 8X U A - 5X 12 - A T X8A - C T 8C = 

A T 8X 22 A - 5X 22 + 8A T XA + A T X8A + C T 5C + 5C T C = 



From (1100b . we have that 

8X n = 

8X 2 2 = —8X12 — 8X^ 2 
Thus for the second term of (|9TT ), we have that 

tr [{B + 8B) T {P + (JPi + 5P 2 + o(8P 2 ))(B + SB) 
=tr [B T PB] + tr [8B T PB + B T 8P 1 B + B T P8B] 

+ tr [8B T P8B + 8B T 8P 1 B + B T 8Pi8B + B T 8P 2 B] + o(5P 2 ) 
It is readily verified that in (11021 ), 
tr [P T PP] = 

tr [8B T PB + B T 8P 1 B + B T P8B] = 
tr [8B T P8B + 8B T 8PiB + B T 8P 1 8B + B T 8P 2 B] 
< \\8B\\ 2 F \\P\\ F + 2||P|| F ||5P|| F || ( 5P 1 || F + ||B|| F ||5P 2 || F 



(99) 



(100) 



(101) 



(102) 



(103) 



Next, from (1931 ) and (19~6T ), we will derive the F-norm bounds of <5Pi and 8P 2 , respectively. From 
(1931 ), we have that 

vec(5Pi) =-[A T ®A T - I] _1 [{A T P ® J)^ + (I ® A T P)] vec(<L4) 

- [A T <g> A T - /] 1 [(/^ ® C< T ) + (C T ® I 2nx )K ny> 2 nx ] vec(8C) 



--Mi 



vec(<5Po) 
vec(8H) 



(104) 
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Thus, with the confidence level a(n 2 y + m 2 n^, xi), we have (18TI ). Likewise, from (|96l ), we have 
that 

vec(5P 2 ) = -(A T ®A T - /)- 1 vec(5A T 5P 1 A + cL4 T P5A + A^dP^A + 5(7 T 5C) (105) 

Taking 2-norm on both sides in (11051) yields ((82l in the deduction of which we use the fact 
1 1 vec ( - ) 1 1 2 = || ■ H-F- From (| 1 02b ~ (11051) and <|9T1 ), we have the % norm model error bound 
given by (f80l) where with the confidence level a(ra 2 + m 2 n^,x 2 ), the F-norm bounds for 5 A, 
5B, 5C and 5F are given by © ~ ®. ■ 

B. Hoo Norm Model Error Bound for the Identified System 

In this subsection, we will propose two approaches, based on perturbation analysis and LMI 
technique, respectively, to computing the 'H CXi norm model error bound for G e (z) with a given 
confidence level. 

Assume the data size N is sufficiently large such that (|24l) approximately holds. Thus, with a 
chosen confidence level a(n 2 e + m 2 nl, Xa), we have the F-norm error bounds of the state space 
matrices in innovations model, shown in (|83l ) ~ (f86l) . 

A straightforward way to derive the "Hoc -norm bound of the error system is by perturbation 
analysis. Asymptotic expansion of the transfer function G e (u) for the original system with respect 
to the parameters for the identified system G e (u) in (l90l) yields 

G e (u) - G e (co) 

={C - 5C) (e iu> I -A + 5A) 1 {B-5B)+F-5F-C (e^I - i)"* B - F 

= - 5C{e iuJ I - A)- 1 ^ - C{e iuJ l - Ay l bA{e iuJ I - A)~ l B 

- C[e iuJ l - A)- l 5B -5F + o(A) (106) 

where A = {\\5A\\ F , \\5B\\ F , \\SC\\ F , \\5F\\ F ) T . Taking Hoc -norm in both sides of (IT061 and 
then applying triangle inequality, submultiplicative inequality (Hoo norm is defined on a closed 
Banach space), and the fact || • || 2 < || • \\ F yields 

\\G e (u) - G? e (w)||Qo 
<\\C(^I - A)- 1 || O0 ||M|| J ,||||(e <w / - A^BWoo + \\5F\\f 

+ \\SC\\ F \\(e^I - Ay'BW^ + \\C(e^I - i^lU^I^ + o(\\A\\ F ) (107) 
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With the identified state-space matrices A, B, C and F and their F-norm error bounds given by 
(1831) ~ (l8~6l) . the ^oo-norm bound of the error system (l90l can be explicitly derived from (11071 ) 
with a given confidence level a(n 2 e + w?v? el x 2 a )- 

Next, we will propose an LMI-based approach to deriving an upper model error bound for 
G e (z). The following Lemma, similar to Lemma 2.1 in [1261 . provides a useful tool to address 
the uncertainty of the state-space matrices in innovations model. 

Lemma 2: Let A £ P nxn , Q = Q T e R nxn , H k £ P nx \ and P fc £ P^' x " (Jfc = 1, • • • , K) be 
given matrices. If there exist K positive scalars //& > (k = 1, • • • , K) and a positive definite 
matrix P > such that 
-P 



PA 




PP 2 - 


• PH K 















-/ii/ 












-//a/ ■ ■ 





(sym) 

















< 



Then the following inequality holds 



K 



(A + J2 H k F k E k ) T P(A + £ P fc P fc P fc ) + Q < 



fc=i 



fc=i 



for each F fe (A; = 1, • • • ,K) satisfying ||Ffe|| 2 < 1. 

Proof: From (I108I ). applying Schur complement yields 



(108) 



(109) 



A T P 



PA 



Q + SfcLi H>kElE k 
and from the right by 

A 

Q + Z)fc=i ^kE\E k 



< 



Multiplying (II 101) from the left by 



p- 1 

/ 



A T 



p- 1 

/ 



< 



Also, for each k, we have that 



Pfc 









(110) 



yields 



(HI) 



>0 (112) 
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H k F k E k 
ElFlEl 



(113) 



which can be simplified as 

HkE k F k F k E k 
The condition \\F k \\ 2 < 1, is equivalent to 

F[F k <I (114) 

From ffTT3T> and (II 141) . in (11111 we have that 

-P- 1 A + J2i =1 H k F k E k 

a t + Ef=i ^ r ^ T ^ r 

which is equivalent to (11091) . ■ 
To compute the Hoc norm of the error system (|90l . we consider the following minimization 
problem 

mhi7 2 
p 



< 



(115) 



s.t. 



A B 


T 


P 




A B 




P 


C T 




I 




C T 




1 2 I _ 



< 



P > 



(116) 



where 



and 



A 



A 




B 










, B = 




, c = 


c -c 




A 




B 









(117) 



< e 3 , II5FIL < e 4 



(118) 



This LMI problem is equivalent to ||G e (o') — G e {uS)^Uoo < 7 f° r an Y ll^ll^ < e i» ll^ll-P < e 2> 
< £3, II^IIf < £4 where ei ~ e 4 are given by the square roots of the right hand sides 
of (1831) ~ (l86l) . Then we will use Lemma [2] to address the uncertainty of A ~ F in (|116l) . For 
brevity, we denote 



(119) 





i 


















, c = 






I 




B 
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and 





















#1 = 










, #3 = 





, if 4 = 
























6 A SB SC SF 

fx — , r 2 — , r 3 — , ^ 4 — 

^2 €3 €4 

E x = y^I 0], E 2 = [0 yfel], E 3 = y^I 0], E 4 = [0 v^J] 
Then the first LMI in (II 161) becomes 



(120) 



( 


A 


B 




C 






fc=i 



)' 


p 


( 


.4 


B 




7 




c 






fe=i 





p 


)- 


^ 2/ . 



<0 (121) 



where it is readily verified that \\F k \\ F < 1 for k = 1, 2, 3, 4. Due to the fact that || • || 2 < || • 
we also have \\F k \\ 2 < 1 for k = 1,2,3,4. By Lemma [2l we arrive at a suboptimal, convex 
minimization problem 



mm 7 



S.t. 



p 




p 




A 


B 




P 


H x ••• 


P 


I 




I 




C 







I 




7 



P4 



p 



7 2 / 







(sym) 






-/i 4 7 



< 



P > 

^ fe >0 (As = 1,2,3,4) 



(122) 



whose minimum 7 2 is an upper bound of the 7 2 in (11161) such that ||G e (w) — G e (u)\\ nao < 7. 

V. NUMERICAL RESULTS 

In this section, we evaluate our theorems by applying them to two numerical examples 
in stochastic subspace identification and H 2 and norm model error bound estimation, 
respectively. 
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A. Example in Stochastic Subspace Identification 

We present Monte Carlo simulation results for the identification of one typical MIMO system 
(i.e., system with its poles very close to the unit circle) where for most of the identifications, 
the DARE © fails. The results shown afterwards are all from the corrected stochastic subspace 
identification in which (fl"5l) ~ (fTTT) are used to impose the positive realness on the covariance 
model. 

The example corresponds to the following MIMO innovations model: 



( : 'k 





0.874 


0.8 




0.18 


0.85 


Xk+l = 






Xk + 








-0.2 


0.96 




-0.25 


-0.4 




-0.3 


-0.65 








Vk = 






Xk + e k 






0.76 


-1.1 









(123) 



0.075 0.037 
0.037 0.068 



and its poles shown in Figure [Q 
1231) for 200 times in each of which the data size is fixed 



with the covariance Q = E[eke^\ 

We identify the simulated system ( 
at VV = 2500 and m in block-Hankel matrix dU) is chosen to be 4. To reflect the performance 
of the proposed identification procedure, in terms of T-L 2 an d Hoc we define two relative errors 
E 2 and Eoc 



Eo 



\G e (e^) - G e (e^)\\ H2 



E n 



\GJe 



GJe 



(124) 



\\Ge{e^)\\u 2 1 \\G e {? 
where G e (e %u> ) is given by (|40l ). Figures |2] shows the two indices E 2 and E^ for 200 identifi- 
cations. In particular, their respective expectations are 0.5055 and 0.6159, and their respective 
standard deviations are 0.1936 and 0.2669. Figure |3] shows the frequency response for one 
identified G e {u) with E 2 = 0.4360 and E^ = 0.5272. The simulation results show that the 
improved stochastic subspace identification procedure works well and guarantees a valid model 
returned even if the system poles are extremely close to unit circle in the case of insufficient 
data. 



B. Example in the comparison of the asymptotic variance of the transfer function 

This example compares the asymptotic variance of the transfer function estimator G e (z) = 
C '(zl — A)~ l B + F ', computed using (|24j) . (l64l) and (|67|) . with the sample variance obtained from 
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Fig. 1. Poles (*) of the innovations model J 1 23t . 



1.5 




Simulation times 



Fig. 2. Ti.2 and T-L^ norm relative error of the identified system for 200 identifications. 
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0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 
normalized frequency co (rad/s) 



Fig. 3. Frequency response for the original G(oj) (solid) and the identified G{ui) (dashed). 

200 Monte Carlo simulations. Consider the following setup: 

A = 0.8, C = 0.1, K = 0.35, Q = 0.001 

, and set the data size to be 10 4 . 

Figure |4]illustrates that the variance computed from the results of this paper agrees well with te 
sample variance obtained from the Monte Carlo simulations. Since the second order statistics is 
used in (|24l) and only the first-order perturbation is used in (l64b and (1671) . the asymptotic variance 
requires more data to obtain an acurrate estimation than its counterpart in the deterministic case 

na. 

C. Example in I-L2 and Hoc Norm Model Error Bound Estimation 

In this subsection, we will first identify a stochastic system using the proposed stochastic 
subspace identification procedure, and then use the estimated state-space matrices to approxi- 
mately compute the I-L2 and Hoo norm model error bounds with a given confidence level. When 
the system poles are very close to the unit circle, according to our simulations a tight model 
error bound with a high confidence level requires an enormous data size that conflicts with the 
practical conditions. Thus, we choose a different simulation example from (1 1 23b . 
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i 1 1 1 1 1 1 1 1 

0.5 1 1.5 2 2.5 3 3.5 

normalized frequency, rad/s 



Fig. 4. Asymptotic variance and sample variance (Monte Carlo estimation) versus normalized frequency (cj G [0, n]): the solid 
line represents asymptotic variance and the dash line represents sample variance. 



Consider an MIMO innovations model: 

Xk+l = 

Vk = 



0.58 


0.23 




0.15 


0.1 






x k + 






-0.39 


0.82 




-0.25 


-0.4 


-0.3 


-0.65 












x k + e k 




0.76 


-1.1 









with the co variance Q = E\e k e^\ 



l-L^ norm ||G e (uO 



Wo 



(125) 



, the V.2 norm ||G e (a;)||^ 2 = 0.5113, and 



0.075 0.037 
0.037 0.068 

0.9774. In the identification procedure, we choose m = 4 in © and to 
get a tight model error bound with a high confidence level, the data size iV is set to be 1 x 10 5 . 
g in (|72|) is a normally distributed vector, i.e. g — > Af(0,I 68 ). We choose a high confidence 



level a(68,88.5) = 95.18%, i.e. with the probability 0.9518, such that g 3 



< 



Xo.9518 



88.5 in 
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( |72l ). By Theorem [3j (11071 ) and (11221) . we compute the H 2 an d "Hoo norm model error bounds 

\\G e {u>) - G e (u)\\ Ha <0.1059 

||G e (w) - G e (u)\\ Hoo <0.0848 

||G e (w) - G e (u)\\ Hoo <0.1467 

respectively, with a confidence level 95.18%, while the true H2 and Hoo norm for the error 
system (f90b is 

||G e (w)-G e (a;)|| Wa =0.0075 

||G e (o;)-G e (a;)|| Koo =0.0164 

It is shown that (11071) is more likely to achieve a tighter Hoo norm model error bound than 
(fl22l . 

VI. Conclusion 

In this paper, a new and straightforward LMI-based optimization approach was proposed 
to impose positive realness on a formerly identified covariance model, guaranteeing a positive 
definite solution to the DARE © and thus a valid innovations model. As can be seen from the 
numerical results, this approach performs well even if the system poles are very close to unit 
circle in the case of insufficient data. 

Later, we presented a complete asymptotic analysis of the stochastic subspace identification 
algorithm. It is shown that, if the data size is sufficiently large, the estimated state- space matrices 
of the covariance model are normally distributed. Also, using Chi-square cumulative distribution 
function, we derived the asymptotic F-norm error bounds of the estimated state-space matrices in 
innovations model. By combining these asymptotic results, the "H 2 an d "Hoo norm error bounds 
for the identified system are explicitly derived with a given confidence level. Numerical example 
in 1-L 2 and model error bound estimation for the identified system is provided to illustrate 
the performance of the proposed estimation approaches. The tools developed in this paper allow 
to extract the weighted additive "Hoc norm model error bound suited for robust controller design 
from system identification, which we postpone to future work. 
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